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This article is a contribution to the understanding of fluctuations in the out of equiUbrium 
dynamics of glassy systems. By extending theoretical ideas based on the assumption that time- 
reparametrization invariance develops asymptotically we deduce the scaling properties of diverse 
high-order correlation functions. We examine these predictions with numerical tests in a standard 
glassy model, the 3d Edwards- Anderson spin-glass, and in a system where time-reparametrization 
invariance is not expected to hold, the 2d ferromagnetic Ising model, both at low temperatures. 
Our results enlighten a qualitative difference between the fluctuation properties of the two models 
and show that scaling properties conform to the time-reparametrization invariance scenario in the 
former but not in the latter. 
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I. INTRODUCTION 

In a series of papers the idea that time-reparametrization invariance (TRI) should be the reason for large non- 
equilibrium fluctuations in glassy systems was introduced. This claim was substantiated with the analysis of the 
Martin-Siggia-Rose-Janssen-deDominicis (MSRJD) action for Langevin stochastic processes [Il-ll and by numerical 
simulations of various glassy models 0, 0, H, [6| . The state of the art of these ideas is summarized in [3] . 

The non-equilibrium situation we are concerned with is a dynamic process where a system is quenched below a 
transition or dynamic crossover with some protocol and it is let relax subsequently. The most interesting time regime 
after the quench is the one in which the size of the sample is the largest (though still finite) scale and times are long 
but not as much as to allow for equilibration. The system is then in an asymptotic non-equilibrium regime in which 
a kind of dynamic scaling is very often realized. In this stage the behavior of time-dependent observables is rather 
well understood 

The sources of fluctuations in the dynamics of these systems are disorder in the initial conditions, thermal noise, 
quenched disorder and external forces, if present. In [7| and references therein the proposal was to average over 
any sort of quenched randomness with the purpose of bringing together disordered systems - such as spin-glasses, 
superconductors, or other - and conventional glasses - of atomic, molecular or other nature. External forces were set 
to zero and had no effect on the dynamics. One is then left with initial condition and thermally induced fluctuations 
only. In order to study their effects, beside the fully averaged value X of physical quantities, one has to consider the 
properties of their fluctuating parts X, to be defined properly, which may be sample dependent. The most complete 
characterization of fluctuations is then obtained through the probability distribution P{X). In particular, since the 
relevance of fully-averaged two-time quantities such as the self-correlation C(tiA2) function, the linear autoresponse 
x(^ij^2) and their fluctuating parts, to the issue of TRI was underlined in [ll,[2|,Q, one is interested in the probability 
distributions of C and x studied as a function of C, ti, t2, the system size and the other parameters in the model. 
This is, however, a formidable task. Less ambitious but still very instructive is to analyze the flrst moments of 
the distributions, such as the generalized variances that we shall define in the body of this paper. This is, indeed, 
the route that we shall follow in this article. Resorting to a Langevin description, in Sec. |TT] we set the formalism 
for continuous variables systems and then present analytic expressions for these variances in terms of the MSRJD 
fields, deriving scaling laws based upon the TRI scenario. Next, recalling previous results on linear and beyond-linear 
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response theory [ll| - [l3| . we relate these variances to quantities that can be measured in Ising spin systems in Sec. IIIII 
In Sec. lIVI we put the TRI predictions to the test in the Edwards- Anderson (EA) model in d = 3 by means of extensive 
numerical simulations. The results, which refine and extend the analysis performed in 0, are consistent with the 
predictions of TRI. In the same section we present another set of numerical simulations performed on a coarsening 
system, the ferromagnetic Ising model in d = 2, for which TRI is not expected to hold [9], showing a different pattern 
of fluctuations, already at a qualitative level. A summary and the conclusions are contained in Sec. |Vl 



II. FORMALISM FOR SYSTEMS WITH CONTINUOUS VARIABLES AND TIME 
REPARAMETRIZATION INVARIANCE 



A. Generating functional 



The stochastic dynamics of continuous variables in a macroscopic system are usually described with a Markovian 
Langevin equation with additive white noise. The ensemble of such coupled stochastic equations can be recast into a 
path-integral called Martin-Siggia-Rose-Jenssen-deDominicis (MSRJD) generating functional that is better adapted 
for a general type of analysis. This formalism has been reviewed in several publications pjO| so we shall not repeat 
the construction of the path-integral here. It reads 

Z['n,ifi] = J Vcj)Vi$ e-^I*'''^'^! Picmx,to)] (1) 

with 

S[<^>;r],rj] ^ / dtdx \ i^{x,t)[(j){x,t) + ^"^^ - h{x,t)] + T{i${x,t)f + ri{x,t)(j){x,t) + fj{x,t)i^{x,t) \ (2) 



S<p{x, t) 

where, to simplify the notation, we focused on a one-dimensional theory described by a scalar field </>. io is the initial 
time, when the system is set in contact with a thermal bath at temperature T. We set fcs = 1 and we absorb the 
friction coefficient into a redefinition of time. Pic [4>{x, to)] is the probability distribution of the field initial conditions. 
Henceforth we shall not write the limits of the time-intervals explicitly. V is the potential energy from which the 
deterministic force in the Langevin equation derives, h is an external perturbation that couples linearly to the field 
(j) in the potential. 77 and 77 are two space and time dependent sources. The path integral runs over all possible 
configurations of the field 4> with initial condition 4){x, to). The measure is defined as Vcf) = Y[k=o 00=1 dtpixajik) and 
we jump over subtle discretization issues that are discussed in many publications by recalling that we do not need 
to worry about any determinant contribution to the action. The integral over the auxiliary field i(j){x, t) is similarly 
discretized, namely T>i^ = O^i 110=1 '^^'^(^a, ^fc)- Without loss of generality for our purpose, we take Pjc to be 
flat. We use a compact notation for the two fields, $ = (0 , icj)), in the action and hereafter we label = and 
$1 — iip its components. Due to the normalization of the thermal noise distribution the path-integral is identical to 
one, 2[0, 0] = 1, in the absence of sources; this feature will be useful later. Noise averages translate into averages 
computed within the path-integral formulation with weight e~^^^'^'^^. 

For definiteness, and to parallel the discussion in Q and references therein (to which we refer for further details), 
let us separate a quadratic term from the non-linear one, Vni, in the potential 

V[4']^^J dxdy Jix, y)4>{x, t)cp{y, t) + F„z [0] (3) 

and assume that the symmetric time-independent exchange matrix J is made of independent elements each of them 
Gaussian-distributed P[J] oc 6^P[~(l/2) / dxdy J'^{x,y)K^^{x,y)] with the connectivity matrix. For a 

lattice model with first nearest-neighbours and no self-interactions, as the 3d Edwards- Anderson spin-glass that we 
studied numerically, K~^{x,y) = J^/z for nearest-neighbours and zero otherwise, with J the exchange constant and 
z the connectivity of the lattice. Quenched disorder is treated on average by exploiting the fact that Z[0,0] — 1. 
Concretely, one averages the dynamic generating functional over the J-distribution with no need to use the replica 
trick. This gives rise to a quartic term in the action of the form 

S'4[4>] = J dxdy J dtidt2 i(j){x,ti) (f>{y,ti) + (f>{x,ti) i(f>{y,ti) K{x,y) i(^(a;, ^2) ^(y, ^2) + ^(a;, ^2) ^2) 

= -IJ dxdy J dhdt2 $'^(a;,ti)$"(2/,ti)if(x,y)$''(a:,t2)$''(j/,i2), (4) 



3 



where a, 6 = 0, 1 and 0=1 and 1 = (here repeated indices are summed). This term can be decoupled by introducing 
four spatially local two-time fields Q"'''{x;ti,t2), a,b = 0, 1, with a Hubbard- Stratonovich transformation. One then 
has 



S'4[$] = / VQcxp 



i J dxdy K-^{x,y) J dtidt2 Q''\x-tiM)Q^"{yM,ti) 
~ j dx j dhdt2 $"(:E,ti)g"^(x;ti,i2)$'(^,<2) 



(5) 



The zero-source generating functional can now be written as Z[0,0] = J VQV^ f,-s[QM ^^gre the action 5[(3,$] 
includes all terms, those in Eq. ([5|) that depend on Q and those in Eq. ^ with V replaced by Vni that depend on <i> 
only. 



B. Relation between composite operators of Q and $ fields 

The TRI arguments which are briefly recalled in Sec. Ill C[ rely on the transformation properties of the fluctuating 
Q fields that, in turn, translate into relations between the average values of products of the Q fields. These objects, 
though, are auxiliary in the calculation and are not directly accessible in a simulation or experiment. In order to 
proceed towards a test of TRI the next step, to which this paragraph is devoted, is to relate the statistical properties 
of the Hubbard-Stratonovich fields Q to those of the field $. In particular, our program is to find relations between 
averages of products of Qs and moments of the $s from which directly measurable correlation and response functions 
will be constructed. To this end we perform an infinitesimal translation of the field Q ^ Q + SQ, which induces a 
change in the action S[Q,^ S[Q,^ + 5S[Q,^], with 



6S[Q,^] ^ - J dxj dhdt2 <^''ix,h)<S>\x,t2)SQ''\x;h,t2) + J dxdy K-^(x,y) j dhdt2 SQ''\x;ti,t2)Q^''{y,t2,ti) 

dxdy K~\x,y) j dhdt2 6Q''\x;ti,t2)5(^'^{y,t2,ti) ■ (6) 

The normalization Z[0, 0] = 1 implies (e^^^) = 1 where (• • •) = J V^VQ ■ ■ ■ e^^^"^'^^. To first and second order in 
SQ one has 



dy K-\x,y) {Q'^^iy-MM)) = {'^%xM)'^\x,t2)) , (7) 



K~\x,u) 5-^5k-c5{ti-U) 5{t2~h) = l^(^''{xM)^\xM)~ j dy K-\x,y) Q^'^iyMM) 



and 



(8) 

respectively. The first equation tells us that the average of the Hubbard-Stratonovich Q field, coarse-grained over the 
region in which K^^{x, y) is different from zero, is equal to the correlations of two $ fields. These in turn are related 
to the usual correlation and linear response functions, as will be recalled in Sec. IIIDl Note that the rhs of Eq. ]7]) is 
symmetric under the exchange of {a,ti) and (6,^2) and this implies that the average in the Ihs should also be. Since 
we are not only interested in the average of Q, but also in their generalized variances (which are given by averages 
of four fields, as discussed in Sec. Ill D[) . the content of the second order equation is also relevant. Indeed, the 
expansion of the product in the right-hand-side of the second equation yields three kinds of terms, one of which is 
precisely the average of the product of four $ fields. Another one is the average of the product of two Qs. The last 
cross term is the average of the product of two $ and one Q fields. Using standard manipulations of the path-integral, 
this term can also be recast into averages over products of Q fields: 

{^''ix,h)^\x,t2)Q'''iv;t3,U)) = {Q'\v-MM)Q''^{x-hM)) -5{ti-h)5{t2-U)5ac5M5[x-v) . (9) 

Putting everything together, Eq. ([S]) allows one to write the average of four <I> fields as sums of terms with averages 
of products of two spatially coarse grained Q fields. 

As it will be discussed in Sec. Ill CI TRI provides predictions on the transformation properties of composite Q-fields 
operators under a certain class of transformations. Because we have just shown above that averages of Q operators 
are equal to moments of $ fields, we arrive at the following important relation: all physical observables which can be 
written in terms of moments of $ fields should inherit the symmetry properties of the Q fields. 
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C. Time-reparametrization 

An approximate treatment of the MSRJD action P, HI, 0] ; inspired by the solution to mean- field glassy models 
and mode-coupling approximations Q, is based on the assumption that it can be separated into two kinds of terms, 
involving fast and slow components of the fields. In the slow contribution time-derivatives are neglected, integrals 
involving products of fast and slow contributions are approximated, and the action is fully written in terms of the 
Q fields only. The a, h indices are linked in such a way that the ensemble of slow terms are invariant under the 
global transformation t — >■ h{t), where h{t) is any positive and monotonic function of time. This transformation is a 
time-reparametrization. In terms of the Q fields it implies 



'^\x- 1, s) ^ Q^^'j^ix; t,s)=( ^) " (^) Q'^\x; t, s) , (10) 



mi 

and, for the averages of Q-field products 

{Q'^^''^{xi;ti, si) ... (3''"^"(a;„;t„, s„)) -J> (Q^'^' (a;i; ti, si) . . . Qta" (a;„; «„)} = 

(Q'^i''i(xi;ii,si)...Q""^"(x„;i„,s„)) . (11) 



dti J \ dsi J \ dtn J \ dSn 

Equation (1111) describes how a generic composite Q-field operator transforms under time-reparametrization. In general, 
quantities characterized by the same set of indexes ai, 6i, . . . , a„, 6„ transform in the same way under the time- 
reparametrization change. A particular case is that of quantities that are left unchanged by the transformation, 
that is to say, TRI-quantities. These can be obtained by integrating a composite operator such as the ones entering 
Eq. (|lip over each time associated to an Ofc = 1 or 6^ = 1 subscript. Indeed, the prefactors are eliminated by a change 
of integration variable t h{t). For the cases n = 1 and n = 2, recalling Eqs. ([7]) and one has that averages 
of products of 2n <I>-fields are TRI after integration over the times pertaining to the $° fields. The autocorrelation 
function C and the associated linear susceptibility x, well as the second moments of their fluctuating parts, belong 
to this class, as it will be further discussed in Sec. HID II These quantities are the subject of the numerical study 
carried out in Sec. IIIII 



D. Variances 



In this Subsection we define the quantities we study in this paper, both analytically and numerically. 



1. Fully averaged quantities and their fluctuating parts 

The connected two-point two-time correlation function is 

C(f , y; s) = ( (0(f , t) - (0(f , t))) (0(y, s) - (^(y, s))) ) (12) 

where (...)= / . . . e"'^. The averaged linear response of any observable is defined as its variation with respect 

to an applied perturbation, in the limit in which the latter vanishes. The most studied linear response is the one of 
the averaged field itself with respect to a perturbation h which couples linearly to in the Hamiltonian (a magnetic 
field, regarding as a spin). A simple calculation shows that, within the MSRJD formalism, it is simply related to 
an average of the product of the fields </) and icfi: 

^T{(t>{x,t)imti)) . (13) 

Note that a factor T has been added in the definition of the linear response to ease the notation here and in the 
following. The integral of this quantity is the dynamic susceptibility 

x'^'''\x,y;t,s) = / dh R^'''Hx,y;t,ti) , (14) 

J s 

which in the following will be more simply denoted as x, omitting the superscript We define a particular second 
order response, as the variation of the averaged composite field (t){x, t)(/){y, t) with respect to two instantaneous 



R'^^'^\x,y-t,ti) 



T 



S{^{x,t)),, 



5h{y,ti 
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perturbations applied at ti and t2 on the same spatial points on which the fields are evaluated. Again, a simple 
calculation allows one to show that within the MSRJD formalism this quantity equals an average of, in this case, four 
fields in the absence of the perturbation: 



5h{x,ti)6h{y,t2) 



= r2 {cf>{x, t)i^{x, h)cj>{y, t)t^{y, t^)) 



(15) 



h=0 



In general we define generalized responses response of a composite operator containing the product of 

n-fields with respect to m magnetic fields. These are functions of n + m spatial positions and n + m times, but to 
keep the notation simple we will omit writing coinciding arguments (for instance we will write i?^^'^-'(x; t, <i) instead 
of R^^'^^x,x;t,t,ti) and R^'^''^\x,y]t,ti,t2) instead of R^^''^\x,y,x,y]t,t,ti,t2)). Quite generally, a variation with 
respect to the field h can be replaced by a response field icj) to compute higher-order response functions and, according 
to the convention introduced above, we include m factors of T in the definition of the generalized response functions. 

The next step is to introduce a couple of local fluctuating quantities, C{x, y; t, s) and x{^, y; t, s), with the property 
{C{x, y; t,s)) = C{x, y; t, s) and {xix, y; t, s)) — x{x, y; t, s), to be interpreted as the fluctuating parts of C and x- Since 
different fluctuating quantities can share the same average, this can be done in different ways. A detailed discussion 
of this issue is exposed in [T]| . For continuous variables described by a Langevin equation, recalling the formalism of 
Sec. ini one is naturally led to consider 



C{x, y; t, s) = (l){x, t)(l){y, s) , 

Xix, y;t,s)^T / dti (j){x, t)i(j){y, ti) , 



(16) 



as being the variables whose MSRJD average yield C and x- Still, x involves the auxiliary fleld i(t> and, in consequence, 
it is not directly accessible in an experiment or a simulation. This problem can be bypassed by studying the probability 
distribution of (C, x) through its moments, since, as we show below, these can be expressed in terms of the fleld (j) alone 
relating them to generalized response functions. Starting from the lowest non-trivial moments, with the definitions 
(1161) one can build the variances 



F^^(x,y;t,s) - {5C{x-t,s)5C{y-t,s)) , 
V^^ix;,y;t,s) = {SC{x-;t,s)Siiy;t,s)) , 
V^^{x,y;t,s) = {5x{x;t,s)Sx{y;t,s)) , 



(17) 
(18) 
(19) 



where for a generic quantity X we have defined the fluctuation SX = X — X ^ X — {X). As anticipated above, 
the variances containing 6x can be expressed in terms of the physical field ((> alone, by expressing them in terms of 
generahzed responses. For instance, the variance (fT9l) can be related to a second-order response, as given by Eq. (fTSi) . 

V^^{x,y;t,s) = ^ dtij^ dt2 [R'^^'^\x-,y;t,t,,t2) - R^'''Hi^;t,h)R^'''\y;t,t2)] . (20) 

Analogously, the co-variance V'^^ can be expressed as 

V^^ix,y;t,s) = T f dti { [^{x,t)(l,{x, s) - {(^{x.t^ix, s))] [^iy,t)i${y,^^^ 

J S 

dti [R'^^^^\x,y;t,s,ti) ~ C{x;t,s)R'^^''\y;t,ti)], (21) 



where C{x;t,s) — C{x,x]t, s) is the autocorrelation function and i?*^'^'^^ (x, y; s, ti) — {(j){x,t)4>{x,s)(j){y^t)i(i){y,ti)). 
In this paper, we enforce these relations to compute the moments (|17l) - (|19|) numerically. Higher order moments are 
exceptionally demanding to compute, and we leave their study to further works. In order to improve the statistics of 
the concrete numerical measurements it will be convenient to compute the double spatial integral: 



V^J{t,s)^L 



- T-d 



d'^2 



d^y y^^(x,y;i,s) , 



(22) 



where V-^^ is a compact notation to denote the variances F*-^*-^, V'~'^ and V^^, simultaneously, d is the spatial 
dimension and L"^ is the volume of the system. A comment is now in order. In spite of the name we gave to it, V^-^ 
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is not the variance of a physical quantity; indeed, being hnked to response functions, it can even take negative values. 
Still, in the following we shall broadly use the term co-variances (or simply variances) when referring to y*^*^, V'-^^ 
and yxx. 

Let us now come back to what TRI implies for the quantities introduced above, following the prescriptions of 
Sec. Ill Cl We start from the usual correlation C{x,y;t, s) and linear response x(a;, i, s). The former is written in 
terms of two fields, in terms of Q's it is then a function of (5"° only and it is invariant under Eq. (fTO|) . x involves, 
besides one field (p, also one field, but there is an integral over the time argument of the latter; thus x does not 
transform either. Proceeding analogously for the variances, V'~^'~^ involves only </> fields, V'-'^ and V^^ involve one and 
two fields respectively, together with one and two integrations over the corresponding time variables. Then, also 
the variances are TRI, similarly to C and x- 



2. An equilibrium sum rule 



A sum rule that generalizes FDT to higher powers is easily proven using the invariance of the action and the measure 
under the transformation (j)(x,t) — (j){x^ —t) and itj){x,t) — iij){x, —t) + f3dt(j>{x, —t) in equilibrium fl4|. Indeed, 







d'^x Sx{x; -t, -s) 



f^x {Sxix; t, s) - 5C{x; t, t) + 5C{x] t, s)j 



(23) 



for t > s. The first identity is due to causality and the second one is derived by using the symmetry. As a special 
case one has 



rCx 



(24) 



valid for n — 2 and for variables such that SC{x;t,t) — 0, e.g. Ising spins. 



E. Restricted averages 

For a generic fluctuating observable X, besides the standard full average X — (X) we also introduce the restricted 
average X^ = (X)^, namely an average taken only over the instances (realizations) with a given value C of the 
fluctuating autocorrelation, as suggested in [l^. Clearly one has 

^ = E^c^(^)' (25) 
c 

where P{C) is the marginal probability that a chosen dynamical trajectory has autocorrelation C. For a co-variance 
between two quantities A and B one has 

= V^^P{C) + V^^ (26) 
c 

where V^^ = {AB)^~ {A)^{B)^ is the restricted covariance of 1 and B (for C flxed), and V^^ = Y.^A^B^P{C) - 
A-;^P(C) J2r ^r-Pi^) is covariance of the restricted averages A^ and B^ as C is varied. 



F. Scaling 



In the late stage of the evolution of slowly relaxing systems, some sort of dynamical scaling is usually obeyed. 
Indeed this is what we will find in our simulations. In this Section, we discuss the implications of TRI in such 
systems. Let us start considering a couple of generic local observables A{t,s), B{t,s), such as the autocorrelation C 
or the autoresponse x, with scaling expressions 



A{t,s) = m'-fA 



m 



(27) 
(28) 
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where ^{t) is a certain function of time. We emphasize that we always refer to the slow degrees of the dynamical 
process and that any explicit fast contribution must be subtracted out if present. In some cases, as in coarsening 
systems, ^ can be readily interpreted as a growing correlation length, although this is not crucial in what follows. 
Clearly, if the two quantities transform in the same way under TRI it must he bA — and /a^z) oc in the 

double limit s — ?> cxd and z — > oo, i.e. t ^ s. 

If one time, say t, can be eliminated in favor of B, the quantity A can also be expressed as A(t^ s) = A{B, s) = 
?(s)'''*/A_B['?(s)^''^i?], where Jab = IaJb^, and one can introduce the slope 



Xab{B,s) 



dA{B,s) 



z=^{s)-''BB 



dB dz 
Under time-reparametrization t — > h{t) one has 

A{t,s) ^ ATR(t,s)^£,[Ks)f^ Sa 
Defining Atr similarly to what we did for A, one has Atr{Btr,s) — ^{h{s))''^ Jab [^{h{s))~'"^ Btr] , and 



mt)) 



Xa {Btr, 



dAxRiBxR, s) 



dB' 



TR 



= m^)) 



-bB+b/ 



dfAB{z) 



dz 



(29) 

(30) 
(31) 

(32) 



Z=e,{h(s))-^BBTR 



Comparing Eqs. p9p and ([32]) we conclude that, if the two quantities transform in the same way, namely &a = ^s, 
the slope of the original parametric curve, xab{B,s), evaluated at a value B is equal to the slope of the trans- 



formed parametric curve, X', 



{Btr,s), evaluated at Btr = [£,{s) / S,{h{s))] "^B. In other words, under time- 



reparametrization B and A are transformed in such a way as to preserve the original curve with slope Xab , although 
different parts of the curve are related by the transformation. Considering the self-correlation and the self-response, 
A = X and B = C, X^c turns out to be the so-called fluctuation-dissipation ratio X = T/Teff, Teff being usually 
denoted as the effective temperature (here we simply quote the term, without entering into the delicate problem of 
its interpretation as a physical temperature) [l^. In this case, since according to the predictions of TRI x and C 
are both invariant under time-reparametrization, Eq. (1321) implies a flnite T^ff. On the other hand, in systems where 
TRI is not obeyed C and x may transform differently and this prevents T^ff from being finite. This leads to the 
conjecture that TRI may not be obeyed in coarsening systems with an infinite Teff, as it was analytically shown in 
a solvable case in For the Ising model this issue will be discussed in Sec. IIVB^ 

Let us consider now non-local quantities such as, for instance V-^"^ {x, y; t, s). Since spatial translational invariance 
as well as isotropy are respected in the mean, they are actually functions of r = \x — y\ only. Assuming scaling, one 
has 



F^^(x,y;t,s)=e(=s)''-- fyx. 
The k = component of such a quantity is then expected to scale as 



V,llit,s)^L-' 



d^x / d^y V''^{x,y;t,s)=a'^)''+'-- /, 



yXY 

k = Q 



m 



(33) 



(34) 



In Sec. IIVI we shall study the behavior of the variances Vj^, hereafter simply referred to as F'^^, numerically. Since 
TRI predicts that all the V^^ should transform in the same way, we expect them to scale as 



XY 



m 
[m 



(35) 



with a unique exponent bxY — b. 



III. SPIN SYSTEMS 



In this Section we adapt the definitions of the variances and restricted averages presented in the previous section 
to discrete variables systems. 
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A. Basic definitions 



1. Models 



We consider a spin system described by the Hamiltonian 

H{a) = -'^J^ja^aJ, (36) 

where ai = ±1 are spin variables located on the N = L'^ sites i of a cubic lattice of linear size L. a denotes a 
spin configuration and Jij are the couplings between nearest neighbor sites (ij). A constant Jij — J defines the 
ferromagnetic Ising model while in the EA model we choose random bimodal couplings, Jij — ± J, with zero mean. In 
the following we shall assume J = 1. The evolution of the system is governed by the transition rates w{a'\(7) for going 
from a configuration a to another a' . In the following we consider the dynamics whereby single spins are updated, 
namely w{a'\a) = (1/-/V) ''J^ii'^'W)- is the transition rate for flipping ai, and a and a' may differ in the i-th spin 
only. We use Wi of the Glauber form Wi ~ (1/2) [1 + tanh(— A£^/2T)], where AE is the energy variation due to the 
proposed spin flip and T is the temperature of the heat bath. 



2. Averaged two-time functions 



The autocorrelation function is defined as 

Cit,s) = {a,{t)a,is)). (37) 



Here we have implicitly assumed {<Ji{t)) = Vt, and (...) = (...),/ denotes, for the spin glass model, the double 
average over thermal histories and initial conditions (denoted by (...)j) and over the quenched realization of the 
random couplings Jij (denoted by '^^); for the ferromagnetic Ising model only the first average is present. The 
associated susceptibility is 

X{t,s) = / dh R^^'^\tM), (38) 



where 

(39) 



dhi{s) 



h=0 



is the impulsive auto-response function and hi{s) is an instantaneous magnetic field switched on and off at time s < t. 
These averaged quantities do not depend on the position i and we thus omit the i subscript. 



3. Fluctuating two-time quantities 

In order to study fluctuating two-time quantities in spin models and to put the TRI to the test, we have to provide 
expressions for the variances (|17p - (jl9p in systems with discrete variables, where the path integral representation 
discussed in Sec. Ill AI is not available. Following this can be done by using an extension of the fluctuation- 
dissipation theorem Jj^ which allows one to relate response functions of any order to correlation functions in systems 
out of equilibrium. To linear order one has 



»(i,.s) 



(Ti{t)ai{t) - ai{t)(Ji{s) - (Ji{t) I dtiBi(ti) 



(40) 



where Bi = — [ai — a'i\'Wi{a'\a). It easy to show that this operator rules the evolution of the local magnetizations, 
namely d{ai{t)) / dt ~ {Bi{t)). It can also be shown [l^ that the expression (|40l) holds also in a Langevin description 
if the analogue of Bi, namely the deterministic force —5V/5{(t){x,t)), is used. This strengthens the relation between 
the discrete variable formalism of this section and the one introduced in Sec. |lll With this fluctuating susceptibility 
and Ci{t, s) — ai{t)ai{s) one readily builds V'~"'~" and V'~"^. The procedure to determine the expression of V^^ is more 
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subtle as we now discuss. First, this quantity is related to the second order response in Eq. ([15]), which for discrete 
spins on a lattice is more properly written as 



V , 1, - 5h,{h)5hj{t2 

by 



(41) 

/i=0 



t rt 



F,f(i,s)= / dh / dh R^,'^\t,h,h) - R^'''Ht,h)R('''\t,h) . (42) 



Proceeding as before, we express R^^'^"> in terms of correlation functions through the extension of the fluctuation- 
dissipation theorem mentioned above. Using the expression for the second order response p2)) derived in one 
arrives at [ll| 

-i^ dti {a,{t)a,{t)B,{h)[aj{t) - a^is)]) + ^ dt, d^ {a,{t)(j,{t)B,{t^)Bj{t2)) 

+ ^ dh {a,{tfB,{h)) - x{t, s)x{t, s), (43) 

where B^ = — [cr' — ai]'^w{a'\a). In this way all the variances are written in terms of the spin configurations and 
they can be computed in a numerical experiment. This method to compute response functions in the unperturbed 
evolution [i.e. Ens. (|40)) and (l43l) ]. besides the advantage of having the vanishing perturbation limit built in, is by far 
more efficient [iTj as compared to the standard method where the perturbation is actually switched on, particularly 
for higher order responses, as discussed in [18]. The accuracy of the results presented in Sec. IIVI would not have been 
possible without resorting to this technique. In the following we focus on the variances' k = component, defined 
analogously to Eq. ([22)) . 

^'''^(^'^) = ^E^#''(^'^)' (44) 

where to ease the notation we omitted the under-script A; = 0. Notice that V^^ is also (apart from an overall factor 
TV) the variance of the space average of the autocorrelations, Cr = (1/-^) J2i C'i, and it is hence bound to be positive. 
This quantity is related to the one measured in Q with the only difference that here we take the coarse-graining 
linear size, £, to be equal to the size of the system, £ — L. We recall that in equilibrium the variances are related by 

Eq. ([MD 

4- Restricted averages 

The definitions in sec. Ill El can be straightforwardly applied to problems with discrete variables. 

IV. NUMERICAL SIMULATIONS 

In this section we present the results obtained with Monte Carlo simulations of the 3d Edwards- Anderson model 
and the 2d ferromagnetic Ising model. We choose to present, first, the bare data for the 3d EA model and their scaling 
analysis. Only later we interpret these forms in terms of TRI. While doing so we confront to the 2d Ising model data. 

A. 3d Edwards-Anderson spin glass 

1. Equilibrium behavior 

In this section we study the ^'''^^s when the spin glass is in equilibrium at T = 0.8. The knowledge of their behavior 
will be useful in the interpretation of the non-equilibrium results, described in Sec. IIV A 21 With this perspective in 
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mind, we selected very small sample sizes to allow for a rapid equilibration. We simultaneously monitored finite-size 
effects as explained below. These conditions led us to work with L < 5. In the following, time is measured in Monte 
Carlo steps (MCs). 

The qualitative behavior of the V^^s and the role of the system size can be understood by inspection of Fig. [T] 
Equilibrium dynamics is stationary so the three depend on the time-difference only. V'~"~' is positive, for the 

reasons discussed in Sec. IIII Al Initially (soon after a microscopic time) it grows as a power-law, V'~"~"' ^ {t — 
s)", and after a characteristic time difference t — s ^ tniL) it flattens and attains a constant value V^*-^(L) = 
(1/A^) ^j^- {{(Ji(jj)jY. This shows that the L-dependence of V^'~' {L) is a finite size effect due to L < ^eg. Because of 
the rapid increase of the characteristic time t^/ with L, the asymptotic value is fully reached only for L = 2 in our 
simulations. For L = 5 the flattening is completely out of reach. From the analysis of this dynamic curve we extract 
the value a — 0.4 for the power-law exponent. 




10000 



le+05 



FIG. 1: (Color online.) Monte Carlo simulations of the 3d EA model. In the main part of the figure the zero wave-vector 
variances V'^^ are plotted against t ~ s in a log-linear scale, for systems of linear sizes L = 2, 3, 5 in equilibrium at T = 0.8. 
See the key for the color and symbol code. Averages are taken over 2.2x10* realizations of the random interactions and the 
initial conditions. In the inset the data for L — 5 are shown in a log- log scale (we use —V'~^^ in order to display a positive 
quantity). The bold orange straight segment represents the power law {t — s)" with a — 0.4. 

The covariance V^^ is negative. This can be understood since the equilibrium parametric plot xi^) is a curve 
with local negative slope. Apart from the sign, this quantity behaves initially as V'-"-' , exhibiting a power-law growth 
with an exponent that is compatible with the one found for V'~^'~^ . Around the quantity V'~'^ reaches a minimum 
and then grows towards an asymptotic value V^^-^ which can be shown to be V^'"- = in cases in which the system 
is not magnetized, {(Ti) — 11]. Indeed, in the case L = 2, where ti\i is sufficiently short to let the system reach the 
asymptotic value, one has a very small final numerical value of V^^. 

Finally, V^-"^ is initially positive and grows proportionally to the other two variances. Around Im it reaches a 
maximum and then decreases to a negative asymptotic value given by — —V^'" ll]. 

In short, the initial time-difference dependence is given by 

V^^{t, s) ^ V^^{t - s) ^{t- s)" for t - s < tM{L) . (45) 

Let us stress that curves with different L collapse in this power-law regime. Hence, for instance, the curves obtained 
for L = 5 are also representative of a system with L > 5 at least up to times t ~ 10^ MCs. This fact will be used in 
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the next section, where a sample with L = 10 will be studied up to t = 10^ MCs. The sum rule (|24|) is satisfied at all 
times within our numerical accuracy. 




0001 1 ' — ' — ' — ' — 

1 10 100 1000 10000 

t-s 



FIG. 2: (Color online.) The zero wave-vector variances V after a quench of the 3d EA model from T —¥ oo to T = 0.8. 
In the main panel the data are plotted against t — s in a log-log scale for different values of the waiting-time s given in the 
key. The three groups of curves are V^^, —V'^^ and V'^^ from bottom to top, respectively. Averages are taken over 6.3x10^ 
realizations of the randomness and initial conditions. The (purple) dotted curves are the equilibrium behavior in the case with 
L = 5. In the inset the quantity V^^ + 2V'^^ + V'~^^ is plotted against t — s, proving that it vanishes at short time differences. 
Thus the equilibrium relation applies in the quasi stationary regime. 



2. Out of equilibrium behavior: variances 

We now turn to the study of V'^^ in out of equilibrium conditions. We performed a quench from T — >■ oo to the 
working temperature T = 0.8 at time t = and we followed the relaxation of a system of size L = 10 up to 10^ 
MCs. We checked that in the range of times accessed in our simulations the system is free from finite size effects and 
remains out of e quil ibrium. This is probably due to the fact that, although some kind of growing length ^ is present 
in the system (sl. Il9l. [20j . in the range of times considered here it is still of the order of a few lattice spacings (say 3 
or 4 at most). The choice of the system size and range of times is dictated by the necessity to tame the colossal noise 
afflicting y-^-*^ in order to have sufficiently clean data. Indeed, increasing the system size or running for longer time 
enhances the noise. Concerning size, this is due to the fact that only the couples ij within a distance of order ^ in 
Eq. (j44p contribute significantly to the sum, as expressed by Eq. (j33p . while the others introduce noisy contributions 
(in spite of having zero average). Regarding time, as ^ increases the number {L/£^)^ of independent regions decreases 
and this also enhances the noise in the course of time. Furthermore, being two-time correlations, the V'^^ decay to 
zero as i — s increases, making the signal/noise ratio worse. The data presented below are obtained by averaging over 
more than 6.3x10^ realizations of thermal histories, initial conditions and coupling constants. 

The behavior of the variances is shown in Fig. [51 One finds V'~'''~' > and V'~'^ < 0, basically for the same reasons 
discussed in Sec. lIV ATI while V^^ is positive. For small time differences t~s ^ t(s) ~ s one finds a quasi-equilibrium 
behavior: the variances depend on the time difference t — s only, they superimpose on the curves found in equilibrium 
(for the largest size L = 5, see the discussion at the end of Sec. lIV AT]) and the relation (|24p is verified (see the inset). 
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In short, one has 

\V^^{t, s)| ~ - s)" with a ~ 0.4 m the stationary regime. (46) 

For larger values oi t — s the system enters the aging regime in which non-equilibrium effects become important. 
One observes an explicit dependence on the two times, and a departure from the equihbrium behavior. In this regime, 
the V'^^ start growing faster than the power-law characteristic of the equilibrium behavior. For still larger values of 
t — s, V'-^'-^ and V'-^^ continue their growth steadily but in a slower logarithmic way (see Fig. O, ~ {Intf^^. The 
exponents hxY were found to be hcc = 0.86 and hcx = 0-52 in the region ln(i/s) > 1. On the other hand, V^^ 
seem to saturate to a constant value (or there could be a maximum after which it would decay, we cannot exclude 
this possibility). The difference of the scaling function with the ones of the other variables remains, however, only 
logarithmic and, in some sense, marginal. Notice that the equilibrium part of the variances is not negligible in the 
aging regime. Instead, increasing as a power-law, it is expected to be asymptotically dominant (although we cannot 
access such a long time behavior in the simulations), since either converge to a constant value {V^^) or grow 

logarithmically {V^^ and V^^). 




FIG. 3: (Color online.) The bare data shown in Fig. [2] are collapsed by plotting s '^y^^ (with a — 0.3) against t/ s in a 
log-linear scale for different values of a (see the key). The three groups of curves are eV^^, -l.W'^'^ and V'^^ from bottom 
to top, respectively. The factors -1.5 and 6 are used to achieve a better presentation of the data since V'^^ > \V^^\ ^ V^^. 



Interestingly, as shown in Fig. |31 all variances obey the scaling form 

V^^{t, s) ~ s"" fxY{t/s) in the aging regime, (47) 

with the same exponent a ~ 0.3, the same dependence on the two times t and s through their ratio in the remaining 
factor, but different scaling functions. This fact has important consequences towards the tests of the TRI scenario as 
we shall discuss in Sect. IIVBI 
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FIG. 4: (Color online.) Restricted averages in the 3d EA model for different values of C indicated in the figure and 

s given in the key. The data are shown in a log-log scale. The three groups of curves have been vertically displaced for 
convenience. The (orange) straight segment is the power law {t — s)" , with a — 0.4. In the upper-left inset the parametric 
plot of Tx versus C is plotted for s = 5 x 10'' MCs. The straight (magenta) line is the equilibrium behavior Tx = 1 — C. In 
the lower-right inset the same curves of the main part of the figure are plotted with the rescaling of Fig. |3l namely by plotting 
g-ayxx^ with a = 0.3, against t/s (in a log-log scale). The dotted (magenta) line is s~'^V^^, for s = 5 x 10^ MCs. 



3. Out of equilibrium behavior: restricted averages 

The restricted average of the susceptibility, V^^, is shown in Fig. |4] for three values of the two-field product, 

C = 0.85, C = 0.75, and C = 0.65. Interestingly enough, for all C the sign of the restricted variance is negative, 
namely the opposite of the sign of the whole V^^ . From the parametric plot of x versus C (included in the upper 
inset) one concludes that the value C — 0.85 belongs to the quasi-equilibrium regime, where Tx ~ 1 — C, whereas 
C = 0.65 is far into the aging regime. The value C — 0.75 is somehow in between. For the choice C — 0.85 the curves 
for different values of s collapse, showing that V^^ depend only on t — s and, moreover, the data are well described 

by the power-law V^"^ ^ (t — s)" with an exponent that is consistent with a = 0.4 as found for the full variance V^^ 

in equilibrium and in the quasi-stationary regime. Therefore for values of C in the quasi-equilibrium regime of the 
parametric plot x vs. C, the restricted averages have an equilibrium behavior (even for times t — s > s where the 
whole V^^ falls out of equilibrium). Moving towards the aging regime, for C — 0.75 and C — 0.65 the data collapse is 
progressively lost, starting from the largest values of (i — s)/s. Interestingly, while time translation invariance is lost, 
the scaling of the data transform into the non-equilibrium one, Eq. (|47p. found for the totally averaged variances, 
with the same exponent a ~ 0.3. This is shown in the lower inset of Fig. HI where we scaled the restricted averages 
by plotting s^'^V^^ against t/s. The collapse is quite good for C — 0.65 and it gets worse for increasing C. The 

reason why the collapse for C — 0.85, although definitely worse than the one in the main panel, is not totally lost is 
due to the fact that the behavior V^'^ ^ [t — s)" can be written in the scaling form (H71) with a in place of a. Since 
the exponents a and a are not too different one would need a much larger range of s in order to resolve the different 
curves in the lower inset of Fig. 01 Notice that for C small enough, when the restricted averages scale as in Eq. (|T7)) 
the scaling function seems to be similar to the one of the whole V^^^ apart from the change in sign. This statement 
is proven by comparing —s~°'V^^ with s~°-V^^ (dotted line in the inset). 
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B. Numerical results and time-reparametrization invariance 

In this Section we re-consider the numerical results of Sec. lIV Al for the spin glass and we discuss whether, how, and 
to which extent they may be interpreted as supporting TRI in the relaxation dynamics of the 3d EA model. We also 
report, in Sec. lIV ITS] the result of simulations of coarsening systems, where TRI is not expected to hold, to signal the 
differences with the dynamics of the disordered model. 



1. 3d Edwards- Anderson model 

A first observation regards the transformation properties of the variances V^^. We have shown in Sec. HI] that if 
TRI holds the variances must all transform in the same way, analogously to what happens for the usual two-time 
functions C and x- According to the discussion of Sec. Ill Fl when scaling holds Eq. ([55]) must be obeyed with a unique 
exponent axY, independently of the choice of XY. Our findings, expressed in Eq. (j47|) . show that scaling holds with 
an algebraic function £,{t), and that all the variances scale with the same exponent a ~ 0.3. This result, therefore, is 
in agreement with the predictions of TRI. Physically, the scaling (|47|) with a single exponent axv = a amounts to 
say that, fixing t/s, the variances of C and x s-re proportional. This means that basically the fiuctuations of x are 
triggered by those of C (or vice versa). As we shall discuss in Sec. lIVB^ the situation is very different in coarsening 
ferromagnets. 

The behavior of the restricted variances may help us better understand this point. As expressed above, if TRI holds 
the non-equilibrium fluctuations of x ^^'^ C a-^e deeply related, and one can modify the former by manipulating the 
latter. Restricted averages are a tool to produce an extreme perturbation of the fluctuations of C, because we impose 
VS'^ = 0. As discussed in Sec. IIV A^ this radically changes the behavior of the fluctuations of Xj since even the sign 

of VX^ is the opposite of the one of V^^. The change in sign and the fact that V^^ also scales according to Eq. (|47p. 

c c 
imply that 

>^y|>^P(C), (48) 
c 

from Eq. (p6| . This means that the most relevant source of fluctuations of x are the fluctuations of C, while the 
fluctuations of x which are independent of the variations of C are less important. Again, this property is not found 
in coarsening systems, as we shall show in Sec. IIVB^ 



2. Comparison with coarsening systems 

TRI is not expected to hold in coarsening systems quenched below Tc Q where an infinite effective temperature [lB| 
is developed. While the features of the 3d EA spin-glass discussed insofar basically agree with what one would expect 
in the presence of TRI, in this Section we show that the same analysis performed on a phase-ordering system, the 
Ising model, leads to a striking disagreement. 

The behavior of the variances in the Ising model quenched below Tc was studied numerically in fll| . A scaling 
form analogous to Eq. (j47p was found, with exponents compatible with acc — d/2, acx — '^xx ~ in d = 1,2 
(and conjectured to hold for higher dimensionality). The scaling functions fxviz) decay as z^^'^ with Acc — 0, 
Xcx = 0; and A^^ = 1/2. This shows that, except at the lower critical dimension d — 1 where the anomaly of a finite 
effective temperature is found [2l|, the variances scale with different exponents. Moreover, the behavior of the scaling 
functions is also very different. As discussed in Sec. [TTl this excludes the possibility of TRI. 

The restricted averages were not studied in phase ordering kinetics so far. In order to establish how generic the 
behavior observed in the 3d EA model is, we computed these quantities in coarsening systems. We simulated the Ising 
model on a c? = 2 square lattice, after a quench from infinite temperature to T = 0. The system size is L = 500. As 
compared to the spin-glass, this relatively large value is needed to avoid finite size effects, since the size of correlated 
regions grows much faster in this problem. As a consequence, large fiuctuations of C around the average are very 
unlikely. Then, when computing V^^ we can collect a significant statistics only in a narrower time interval than in 

the spin-glass case. This can be seen in Fig. [SJ where the behavior of V^^ and V^^ are shown. For V^^, we basically 

reproduce the results in [Tlj (where however a quench to T = 1.5 < Tc was considered). This quantity is positive 
and behaves as V^^ ~ s^/^/^^(t/s), where fxxi^) a function which grows as z^^'^ for large z. Interestingly, the 
restricted average VX^ not only has the same sign as the global average, but it superimposes almost exactly on V^^, 
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FIG. 5: (Color online.) Fluctuating dynamics in the 2d ferromagnetic Ising model. The restricted averages V^^ for different 

values of C as a function of t — s with fixed s — 10^ MCs are shown with (colored) line-points following the code given in the 
key. Similar results are found for other values of s. The data are shown in a log-log scale. The fully averaged variance V"^^ is 
displayed with a solid black line. Averages are taken over 1.1 x 10^ realizations of the randomness and initial conditions. The 
aging regime starts at C ~ 0.65 where V^-^ ~ V^^. 



for t — s sufficiently large (in the aging regime). This means that in the non-equilibrium regime V^^ is completely 
determined by V^^. According to Eq. ([26l) . this implies 



(49) 



which should be compared to the property (|48|) found in the 3d EA model. This indicates that the non-equilibrium 
fluctuations of x are entirely produced at fixed C, and are then independent of those of C. This behavior, which is 
qualitatively different from the one observed in the 3d EA model, is a clear demonstration of the lack of TRI in this 
system. 



V. SUMMARY AND CONCLUSIONS 



This paper is devoted to the study of the properties of fluctuations in aging systems. Specifically, we concentrate on 
the fluctuations of the autocorrelation and dynamic susceptibility, since the relevance of their relations to the issue of 
TRI was pointed out in [l|,|3,0. In the MSRJD path-integral formulation of Langevin processes, which is overviewed in 
the first part of the paper, response functions are associated to the auxiliary field icj) and hence their fluctuating parts 
cannot be expressed in terms of the physical field alone. A possible way out of this limitation is to consider moments 
of such fluctuating parts, since these can be related to generalized response functions which are directly measurable. 
In particular, in this paper we focused on the variances V'^^ , because higher moments are numerically out of reach 
even with the highly efficient method adopted here, based on the generalization of the fluctuation-dissipation theorem 
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out of equilibrium and beyond linear order. The behavior of these variances fits in the general framework of two-time 
quantities in aging systems. For short time differences a quasi-equilibrium stationary behavior is observed, while in 
the aging regime time-translation invariance is lost and is substituted by a different, scaling symmetry with a unique 
exponent axY — o, for all the variances. In the analytic part of the paper we derived the scaling properties of the 
variances and we showed that such a unique exponent is a prediction of TRI. The numerical results in the Sd EA 
therefore gives support to TRI. 

In order to better characterize the fluctuations, we also computed the restricted average V^'^- The values of the 

restriction C is an additional parameter that tunes the crossover from the quasi-equilibrium behavior to the fully 
aging regime. In particular, the latter is obtained letting C be small (with respect to the EA order parameter). In 
this case, the restricted average scales as the fully averaged one, but its actual value is not the main contribution to 
the latter. These results have to be compared to the ones found in the Ising model. In this clean coarsening case the 
exponents axY do depend on the fluctuating observables and the restricted average V^'^ is of the same order as the 

full average t^-*^-^ , showing that fluctuations of C are not the main source for those of at variance to what it is found 
in the disordered model. These qualitative and quantitative differences between these two systems add to differences 
found in the relation between the averaged linear response and the self-correlation, i.e. the effective temperature jj^. 

Our results suggest that the fluctuating out of equilibrium dynamics of the 3d EA model and of coarsening systems 
are different and that this can be related to the fact that TRI is realized in the former model while it is not in the 
latter. This should also be linked to the fact that the effective temperature is finite in the former and infinite in the 
latter cases. 

These ideas should be checked numerically and analytically in other out of equilibrium systems. Natural cases to 
analyze are the dynamics of atomic glass models (like soft sphere, Lennard- Jones mixtures or other) with molecular 
dynamics or the evolution of kinetically constrained models [2^ . The Glauber Ising chain is an example that could 
be dealt with, possibly, analytically, exploiting the techniques derived in [2^. Studies of critical dynamics '24] and of 
the relaxation of a model at its lower critical dimension should help to complete this picture. 
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